home *** CD-ROM | disk | FTP | other *** search
/ The X-Philes (2nd Revision) / The X-Philes Number 1 (1995).iso / xphiles / hp48hor1 / poly.src < prev    next >
Text File  |  1990-10-20  |  6KB  |  220 lines

  1. %%HP: T(3)A(D)F(.);
  2. DIR
  3.  
  4. @ POLY, a polynomial utilities collection.
  5. @ by Wayne H Scott
  6.  
  7.   ZTRIM  @ remove leading zeros
  8.     \<< LIST\-> \-> n
  9.       \<< n
  10.         WHILE ROLL DUP 0 ==
  11.         REPEAT DROP n 1 - DUP 'n' STO
  12.         END n ROLLD
  13.         IF n 0 ==
  14.         THEN { 0 }
  15.         ELSE n \->LIST
  16.         END
  17.       \>>
  18.     \>>
  19.  
  20.   RDER  @ derivative of a rational function
  21.     \<< \-> F G
  22.       \<< G F PDER PMUL G PDER { -1 } PMUL F PMUL PADD G G PMUL
  23.       \>>
  24.     \>>
  25.  
  26.   IRT  @ invert root program (roots -> polynomial)
  27.     \<< LIST\-> \-> n
  28.       \<<
  29.         IF n 0 >
  30.         THEN 1 n
  31.           START n ROLL { 1 } SWAP NEG +
  32.           NEXT
  33.         ELSE { 1 }
  34.         END
  35.         IF n 1 >
  36.         THEN 2 n
  37.           START PMUL
  38.           NEXT
  39.         END
  40.       \>>
  41.     \>>
  42.  
  43.   PDER  @ derivative of a polynomial
  44.     \<< LIST\-> \-> n
  45.       \<< 1 n
  46.         FOR i n ROLL n i - *
  47.         NEXT DROP
  48.         IF n 1 ==
  49.         THEN { 0 }
  50.         ELSE n 1 - \->LIST
  51.         END
  52.       \>>
  53.     \>>
  54.  
  55.   PF  @ partial fractions; handles multiple roots
  56.     \<< MAXR \->NUM { } \-> Z P OLD LAST
  57.       \<< 1 P SIZE
  58.         FOR I P I GET \-> p1
  59.           \<<
  60.             IF p1 OLD \=/
  61.             THEN Z p1 EVPLY 1 P SIZE
  62.               FOR J
  63.                 IF P J GET P I GET \=/
  64.                 THEN p1 P J GET - /
  65.                 END
  66.               NEXT p1 'OLD' STO { } 'LAST' STO
  67.             ELSE
  68.               IF { } LAST SAME
  69.               THEN 1 { } 1 P SIZE
  70.                 FOR K P K GET
  71.                   IF DUP p1 ==
  72.                   THEN DROP
  73.                   ELSE +
  74.                   END
  75.                 NEXT IRT Z SWAP
  76.               ELSE LAST LIST\-> DROP
  77.               END RDER DUP2 5 PICK 1 + 3 ROLLD 3 \->LIST 'LAST' STO p1
  78.               EVPLY SWAP p1 EVPLY SWAP / SWAP FACT /
  79.             END
  80.           \>>
  81.         NEXT P SIZE \->LIST
  82.       \>>
  83.     \>>
  84.  
  85.   FCTP  @ factor polynomial
  86.     \<<
  87.       IF DUP SIZE 3 >
  88.       THEN DUP BAIRS SWAP OVER PDIV DROP FCTP
  89.       END
  90.     \>>
  91.  
  92.   RT  @ find roots of any polynomial
  93.     \<< ZTRIM DUP SIZE \-> n
  94.       \<<
  95.         IF n 3 >
  96.         THEN DUP BAIRS SWAP OVER PDIV DROP \-> A B
  97.           \<< A RT B RT 
  98.           \>>
  99.         ELSE
  100.           IF n 2 >
  101.           THEN QUD 
  102.           ELSE LIST\-> DROP NEG SWAP /
  103.           END
  104.         END 
  105.       \>>
  106.     \>>
  107.  
  108.   L2A  @ LIST->ARRY and ARRY->LIST
  109.     \<<
  110.       IF DUP TYPE 5 ==
  111.       THEN LIST\-> \->ARRY
  112.       ELSE ARRY\-> 1 GET \->LIST
  113.       END
  114.     \>>
  115.  
  116.   PADD  @ add two polynomials
  117.     \<< DUP2 SIZE SWAP SIZE \-> A B nB nA
  118.       \<< A L2A B L2A
  119.         IF nA nB <
  120.         THEN SWAP
  121.         END
  122.         IF nA nB \=/
  123.         THEN 1 nA nB - ABS
  124.           START 0
  125.           NEXT
  126.         END nA nB - ABS 1 + ROLL ARRY\-> 1 GET nA nB - ABS + \->ARRY + L2A
  127.       \>>
  128.     \>>
  129.  
  130.   PMUL  @ multiply two polynomials by doing a convolution
  131.     \<< DUP2 SIZE SWAP SIZE \-> A B nB nA
  132.       \<< { }
  133.         IF nB 1 >
  134.         THEN 2 nB
  135.           START 0 +
  136.           NEXT
  137.         END DUP A + SWAP + 'A' STO  @ pad A with zeros
  138.         A LIST\-> \->ARRY B LIST\-> DROP
  139.         IF nB 1 >
  140.         THEN 2 nB
  141.           FOR J J ROLL  @ reverse B
  142.           NEXT
  143.         END
  144.         IF 3 nA nB + \<=
  145.         THEN 3 nA nB +
  146.           START 0
  147.           NEXT
  148.         END nA nB 1 - 2 * + \->ARRY  @ and pad with zeros
  149.         2 nA nB +
  150.         START DUP2 DOT  @ do dot product and then rotate B
  151.           3 ROLLD ARRY\-> SWAP nA nB 1 - 2 * + 1 + ROLLD \->ARRY
  152.         NEXT DROP2 nA nB + 1 - \->LIST  @ put dot products into list
  153.       \>>
  154.     \>>
  155.  
  156.   PDIV  @ divide two polynomials
  157.     \<< DUP SIZE 3 ROLLD LIST\-> \->ARRY SWAP LIST\-> \->ARRY \-> c b a
  158.       \<< a b
  159.         WHILE OVER SIZE 1 GET c \>=
  160.         REPEAT DIVV
  161.         END DROP \-> d
  162.         \<< a SIZE 1 GET c 1 - - \->LIST d ARRY\-> LIST\-> DROP \->LIST
  163.         \>>
  164.       \>>
  165.     \>>
  166.  
  167.   EVPLY  @ evaluate a polynomial at a point
  168.     \<< OVER
  169.       IF DUP TYPE 5 ==
  170.       THEN SIZE
  171.       ELSE SIZE 1 GET
  172.       END \-> a r n
  173.       \<< a 1 GET
  174.         IF n 1 >
  175.         THEN 2 n
  176.           FOR i r * a i GET +
  177.           NEXT
  178.         END
  179.       \>>
  180.     \>>
  181.  
  182.   COEF  @ equation -> polynomial list  (coefficient extractor)
  183.     \<< \-> E n
  184.       \<< 0 n
  185.         FOR I 0 'X' STO E EVAL 'X' PURGE E 'X' \.d 'E' STO I FACT /
  186.         NEXT 2 n 1 +
  187.         FOR I I ROLL
  188.         NEXT n 1 + \->LIST
  189.       \>>
  190.     \>>
  191.  
  192.   DIVV  @ subroutine used by PDIV
  193.     \<< DUP 1 GET \-> a b c
  194.       \<< a 1 GET c / DUP b * a SIZE RDM a SWAP - ARRY\-> 1 GETI 1 - PUT
  195.          \->ARRY SWAP DROP b
  196.       \>>
  197.     \>>
  198.  
  199.   QUD  @ quadratic equation
  200.     \<< LIST\-> \->ARRY DUP 1 GET / ARRY\-> DROP ROT DROP SWAP 2 / NEG DUP SQ
  201.       ROT - \v/ DUP2 + 3 ROLLD -
  202.     \>>
  203.  
  204.   BAIRS  @ Bairstow's method for quadratic factors
  205.     \<< LIST\-> 1 1 \-> n R S
  206.       \<<
  207.         DO 0 n 1 + PICK 0 0 0 4 PICK 5 n + 7
  208.           FOR J J PICK R 7 PICK * + S 8 PICK * + 7 ROLL DROP DUP 6 ROLLD R
  209.             3 PICK * + S 4 PICK * + 5 ROLL DROP -1 
  210.           STEP 3 PICK SQ 3 PICK 6 PICK * -
  211.           IF DUP 0 ==
  212.           THEN DROP 1 1
  213.           ELSE 6 PICK 6 PICK * 5 PICK 9 PICK * - OVER / 4 PICK 9 PICK * 8 PICK
  214.             7 PICK * - ROT /
  215.           END DUP S + 'S' STO SWAP DUP R + 'R' STO
  216.         UNTIL R\->C ABS .000000001 < 7 ROLLD 6 DROPN
  217.         END n DROPN 1 R NEG S NEG 3 \->LIST
  218.       \>>
  219.     \>>
  220. END